Temporal fluctuations of waves in weakly nonlinear disordered media 1 
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Abstract 

We consider the multiple scattering of a scalar wave in a 
disordered medium with a weak nonlinearity of Kerr type. 
The perturbation theory, developed to calculate the tempo- 
ral autocorrelation function of scattered wave, fails at short 
correlation times. A self-consistent calculation shows that 
for nonlinearities exceeding a certain threshold value, the 
multiple-scattering speckle pattern becomes unstable and ex- 
hibits spontaneous fluctuations even in the absence of scat- 
terer motion. The instability is due to a distributed feed- 
back in the system "coherent wave + nonlinear disordered 
medium". The feedback is provided by the multiple scatter- 
ing. The development of instability is independent of the sign 
of nonlinearity. 



I. INTRODUCTION 

Scattering of waves in disordered media has proved 
to be a nontrivial topic possessing intriguing and still 
not completely understood features JjT3 Accordingly to 
the strength of disorder, one observes a variety of propa- 
gation regimes ranging from ballistic transport, through 
single scattering and wave diffusion, to the Anderson lo- 
calization. In this paper we are interested in the regime 
of wave diffusion, corresponding to a relatively strong dis- 
order, which is, however, still largely insufficient to bring 
the system to the localization transition (kl 3> 1, where 
k is a wave number in the medium, and £ is a mean free 
path) . 

It is well-known, that multiple scattering of coher- 
ent wave in a disordered medium results in a compli- 
cated spatial intensity distribution J(r, t) known as a 
"speckle pattern" . The speckle pattern is highly ir- 
regular and appears random to the eye. It is now 
well established that -the speckle pattern exhibits large 
intensity fluctuationsclu (SI(r, i) 2 ) ~ (7(r,i)) 2 , origi- 
nating from the interference of partial waves arriving 
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at r with completely randomized phases. Here the 
angular brackets (• ■ ■) denote ensemble averaging, and 
6I(r,t) — I{r,t) — (I(r,t)). Besides, the speckle pat- 
tern possesses nontrivial long-range spatial correlation 
Q/lr, Ar) = (SI(r-Ar/2,t)5I(r + Ar/2,t)} even for 
Ar > £_ Xhis correlation is due to interaction of diffusing 
modesJJli 2 ] If the points r ± Ar/2 are far enough from 
the boundaries of the medium, Cgj oc 1/Ar. Recently, 
it has been shown that in a particular case of a point 
source of waves embedded inside a disordered medium, 
there exists an infinite-range contribution to Cgi (r, Ar) 
originating from scattering events taking place in the im- 
mediate neighborhood of the source £3 This contribution 
is highly sensitive to the short-distance properties of dis- 
order, as well as to the source size and shape.li-3 

If the scatterers in the medium are allowed to move, 
7(r, t) fluctuates with time, and the statistics of these 
fluctuations is also a subject of active research. A 
wave propagating in a disordered, multiple-scattering 
medium undergoes a large number of scattering events, 
and hence the scattered intensity is highly sensitive to 
displacements of scatterers. liailH Consequently, the de- 
cay of the intensity autocorrelation function Csi{t,t) — 
(6I(r, t)SI(r, t + t)) is— considerably faster than in the 
single-scattering caseJiatl] Recently, long-range auto- 
correlation-function of intensity fluctuations has been 
measured,Ej and the existence of the universal conduc- 
tance fluctuations (analogous to that in disordered con- 
ductors) has been demonstrated^ 2 ] for optical waves. 
Theoretical analysis of the temporal correlation function 
of multiple-scattered wawes has been extended to ampli- 
fying disordered media,tJ as well as to the case of in- 
tense incident waves producing flows of scatterers in the 
disordered medium.cj'Ell An additional contribution to 
Csi(t, t), originating from scattering in the immediate 
neighborhood of source and/or detector, decaying much 
slower than all the. previously known contributions, is 
predicted to exist JiJ 

High sensitivity of multiple-scattering speckle pat- 
terns to scatterer motion gave rise to a new technique 
for studying the scatterer dynamics in disordered, tur- 
bid media—tlie so-called "diffusing-wave spectroscopy" 
(DWS).Ej'E3l3 The latter is now.— widely applied, in. 
concentrated ,-eolloidal pj.spcnsions,li2lE^rE°l jJpajps,ElTL3 
emulsions granulaiEjeil and biologicalEja media. 
Besides, the DWS has been extended to macroscopi- 
cally heterogeneous turbid media, providing a tool for 



1 



imaging of dynamic l^etogeneitiesETED and visualiza- 
tion of scatterer flowsCJu_j in the bulk of the medium. 
A generalization of DWS technique has been_also accom- 
plished for anisotropic disordered media.KlJTjj Recently, 
the DWS approach has been extended to nonergodic tur- 
bid media-EK 

The above-mentioned, extensive studies of temporal 
fluctuations of multiple-scattered waves, as well as the 
numerous application of DWS, are all restricted to lin- 
ear disordered media. In general, little information is 
available on the subject of multiple scattering in non- 
linear disordered media. Meanwhile, the question con- 
cerning the way in which the nonlinearity affects the 
multiple-scattering speckle pattern still remains open and 
continue to attract the researches. Considerable efforts 
have been made to understand the phenomena of coher- 
ent backscatteriag in disordered media with Kerr-type 
nonlinearity:E3~Ell a narrow dip has been predicted to 
appear on the top of the backscattering peak. Weak lo- 
calization effects are shown to exist in th£_Ladiation of 
second harmonic and difference frequpjcyjEj'ETEj though 
their experimental observation failed. EU Also studiedpwith. 
account for disorder is the optical phase conjugation.EHTE3 
More recently, correlations in transmission and reflec- 
tion coefficients of second harmonic waves have beep in- 
vestigated both theoretically and experimentally,^ and 
the angular correlation functions of fundamental wave in 
a disordered medium with Kerr-type nonlinearity have 
been calculated. O Despite the fact that theoretical de- 
scription of wave scattering in nonlinear media is com- 
plicated by the simultaneous presence of both disorder 
and nonlinearity, the standard diagram technique for im- 
purity scattering has been extended to the case. of. disor- 
dered medium with nonlinearity of Kerr type.BS 

Very recently, it has been shown that the speckle pat- 
tern resulting from the multiple scattering of coherent 
wave in a nonlinear disordered medium with Kerr-type 
nonlinearity, should he extremely sensitive to changes of 
scattering potential,L2l i.e. much more sensitive than the 
linear speckle pattern. This high sensitivity has been ex- 
plained bi2-,the multiplicity of solutions of nonlinear wave 
equation.LJ The multiplicity of solutions has been then 
shown to lead to the temporal instability of the multiple- 
scattering speckle pattern in nonlinear medium, result- 
ing yo-spontaneous fluctuations of scattered wave with 
time.O An important prediction of Ref. [n] is that the 
nonlinearity should exceed some threshold value for the 
instability to develop. The threshold value is principally 
determined by the absorption length L ai or by the sam- 
ple size L, if L < L a . The striking feature of the es- 
tablished result is that the threshold value of nonlinear- 
ity tends to zero in an unbounded medium without ab- 
sorption. Purely elastic, unbounded nonlinear multiple- 
scattering systems are therefore always unstable. The, 
physical origin of the instability is easy to understand.Lil 
Nonlinearity modifies the phases of partial waves propa- 
gating in the medium. The phase modifications are pro- 
portional to the intensity I(r,i) and affect the mutual 



interference of partial waves. As it is this interference 
which is responsible for I(r,t), a sort of feedback estab- 
lishes in the medium. A small modification of J(r, t) 
causes modifications of phases of partial waves which, in 
their turn, produce changes of I(r,t), and so on. It is 
well-known that nonlinear wave systems with sufficiently 
strong, positive feedback may become unstable.EaO As 
an example, we cite a family of nonlinear optical systems 
with two-dimensional feedback,Oil3 where spontaneous 
formation of complicated spatial structures is observed. 
Despite the absence of disorder, such systems exhibit 
transition to seemingly chaotic dynamics with increas- 
ing nonlinearity I 75 " 76 ! An analogy can be drawn between 
the nonlinear optical systems with two-dimensional feed- 
back and nonlinear disordered media by considering the 
scattering as a (three-dimensional) feedback mechanism. 
In the case of disordered media, however, the feedback is 
of random nature and it is therefore hopeless to expect 
regular spatial structures to form. Meanwhile, the insta- 
bility can manifest itself in spontaneous fluctuations of 
speckle pattern. In order to clarify the issue of instabil- 
ity of speckle patterns in nonlinear disordered media, we 
consider the following questions: 

(a). Can the multiple scattering provide a positive feed- 
back mechanism for waves propagating in a nonlin- 
ear disordered medium? 

(6). If "y es '\ how strong the nonlinearity should be for 
the instability to develop? 

A general announcement of our principal answers to-jthe 
above questions has been given in our recent Letter .0 In 
the present paper, we discuss and justify the assumptions 
and approximations underlying our conclusions, provide 
the missing details of calculations, and give a comprehen- 
sive discussion of results. Also developed and discussed 
is the perturbation approach to the calculation of the 
temporal autocorrelation function of multiple-scattered 
wave in a nonlinear disordered medium. It is important 
that the validity condition of the perturbation theory co- 
incides with the condition for the instability threshold 
as obtained by using the self-consistent approach. In 
addition, we give a detailed consideration to an exper- 
imentally important case of moving scatterers, when the 
decrease of the time autocorrelation function is due to 
a combined effect of spontaneous and scatterer-motion- 
induced fluctuations of the speckle pattern. 

The remainder of the paper is arranged as follows. In 
Sec. H we introduce the nonlinear wave equation, and 
discuss how the path integral approach can be applied 
for its analysis. We also formulate the basic models and 
approximations used throughout the paper. Section III 
is devoted to linear disordered media. In this section, we 
provide the expressions for the spatio-temporal intensity 
correlation functions. Although correlations of multiple- 
scattered waves in linear media are well studied nowa- 
days, we present their first, to our knowledge, treatment 
with a simultaneous account for absorption, boundary 
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conditions at the sample surface, and scatterer motion. 
The results of Sec. [II serve as a base for further calcula- 



tions. In Sec. IV, we present a calculation of dephasing of 
waves in a nonlinear disordered medium. Our calculation 
takes into account the fluctuations of the local refractive 
index due to nonlinear effects, as well as the long-range 
spatial correlation of this fluctuations. Three "nonlinear" 
contributions to the dephasing are identified in addition 
to the usual, "linear" term originating directly from the 
motion of scatterers. Further, in Sec. [v| we develop a 
perturbation theory for calculation of the temporal au- 
tocorrelation function of a multiple-scattered wave, and 
show its failure at short correlation times, for sufficiently 
weak absorption. A condition of validity of the pertur- 
bation theory is established by comparing the linear and 
nonlinear contributions to the dephasing found in Sec. 
IV. Section VI presents an alternative, self-consistent ap- 
proach to the calculation of the temporal autocorrelation 
of scattered wave. Development of self-consistent theory 
requires some additional assum ption s which are also dis- 
cussed in this section. In Sec. VII, the main results of 
our self-consistent approach are presented and discussed. 
The multiple-scattering speckle pattern is shown to ex- 
hibit spontaneous fluctuations even in the absence of 
scatterer motion, which we interpret as a signature of 
its instability. A comparison of self-consistent and per- 
turbative results is given, and the condition of the speckle 
pattern instability is shown to coincide with the condition 
of validity of the perturbation theor y. Finally, concluding 
remarks are presented in Sec. VIII. In order to maintain 
the text of the paper readable, we have chosen to collect 
the technical details of calculations in three appendices. 
Appendix [A] is devoted to the derivation of the field-field 
spatio-temporal correlation function. In Appendix RJ we 
compute the spatio-temporal long-range intensity corre- 
lation function. Appendix |c] provides the details of cal- 
culations of path distributions p s (r) and p s (r,r') defined 
in Sec. ffvl. 



II. WAVE EQUATION AND PATH INTEGRALS 

We consider a scalar monochromatic wave of frequency 
u> propagating in a random medium with Kerr-type non- 
linearity. The W&xf amplitude ip(r,t) obeys a nonlinear 
wave equation:OiJj 



■fe(r,i)+e 2 |V>(r,t)| 2 ]}v(M) 



0. 



(1) 



Here fco is the free-space wave number, Eq = Eg-Meg is the 
average (complex) dielectric function, Se(r, t) is the fluc- 
tuating part atthe dielectric function, £2 is the nonlinear 
susceptibilityEj (the two latter quantities are assumed to 
be real). Eq. ([!]) is valid only if fe(r, t) + £2 |^>(r,i)| 
does not change significantly on the time scale of cj _1 . 
The expression in the square brackets of Eq. (pi) can be 



considered as some "effective" dielectric function of the 
medium. General analysis of Eq. (|l|) for arbitrary rela- 
tion between various terms comprising this function con- 
stitutes a formidable task, and is not a purpose of this 
paper. We assume the following hierarchy: 



e 2 |^(r, t)\ 2 « fc(r, t), £ 2 |V(r, t)\ 2 « 4 



(2) 



In other words, we assume that the role of nonlinearity 
is less significant than that of disorder, and that absorp- 
tion is weak allowing multiple scattering of waves in the 
medium. It is then convenient to define the effective re- 
fractive index no = (Eg) 1 / 2 , the absorption length £ a = 
no/(ko£'o), and the nonlinear coefficient n 2 = £2/(2n ), 
which determines the nonlinear correction to the (linear) 
refractive index of the medium: ro(r, t) = hq + n2l(r, i), 
where I{r,t) = \ip(r,t)\ 2 is the wave intensity. 

In this paper, we study the fluctuations of the so- 
lution of Eq. (|I|) with time t. In a linear medium 
(£2 = 0), these fluctuations can only be due to random 
fluctuations of Se(r,t) with time. The fluctuations of 
"0(r, t) are commonly characterized by the autocorrela- 
tion function C^(r,r) = (ifj(r,t)ijj* (r,i + t)) . We as- 
sume that this autocorrelation function is independent 
of t which implies that for fixed r, ip(r, t) represents 
a stationary random process (this is obviously true if 
the sample geometry and the source distribution do not 
change with time, and if Se(r, t) is a stationary ran- 
dom process). Me take Se(r,t) to be a Gaussian ran- 
dom field with zero mean and the correlation function 
C 5e (Ar,T) = (5e(r- Ar/2, t)fc(r + Ar/2, t + r)). For 
a medium composed of point-like scatterers undergoing 
Brownian motion with a diffusion coefficient £>sJl3'EHl 



C 5e (Ar,r) 



4vr/(fc 4 f) 



exp 



AD b t 



(3) 



where k = koriQ, and the mean free path I <C l a is intro- 
duced (a weak scattering limit k£ 3> 1 is assumed). A 
natural time scale for scattering of waves in the medium 
described by Eq. (^) is set by the characteristic time 
needed for a scatterer to move a distance of order of the 
wavelength: To = (4/c 2 Z?s) -1 . From here on, we will be 
interested in short correlation times t <C tq. 

In the linear case (£2 = 0), several approaches have 
been elaborated to analyze Eq. ([j]) . We mention the dia=. 
grammatic techniques, theory of-radiative transfer ,0 
and the method of path integrals E3e3 The three ap- 
proaches are known to give equivalent results for at 
t t . In the present paper, we adopt the method 
of path integrals that was originally proposed in the 
framework of quantum electrodynamics £3 but later has. 
been successfully used in various areas of physics ,E3 
and, in ua*ticular, for the analysis of wave scattering 
problems.E3E3 The method is based on the fact that the 
solution t/j(r,i) of the wave equation (l|) can be written 
in a form of a functional integral, with integration per- 
formed over all possible trajectories (paths) going from 
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the source of waves to r.E 2 ] Since in the weak scattering 
limit {kl 3> 1) different trajectories can be considered 
independent, it appears that the jcorxelation function 
reduces to the following integral:LJ~E3 



poo 

(r, r) = J / P(r, s) exp 
Jo 



da, (4) 



where Iq is the average intensity in a nonabsorbing 
medium, P(r, s) is a weight coefficient of paths of length 
s, and (Atp 2 (r)) s denotes the squared phase difference 
Atp(t,r) = ip(t + t) — ip(t), averaged over various real- 
izations of disorder, and over all possible paths of the 
same length s. From here on, we denote such an averag- 
ing by (• ■ -) a . Note that (Aip(t)) s = forpthe-model of 
Brownian point-like scatterers. Meanwhile JH3cj 



(5) 



where the superscript (0) denotes the linear case. It is 

worth noting that (|A(^ 2 (t))^ ' ) does not depend neither 
on the sample geometry, nor on the source and detector 
positions. Its value is only determined by the scatterer 
dynamics (through the single-scattering correlation time 
to), and the path length s. In contrast, P(r, s) can only 
be calculated if the sample geometry, source distribu- 
tion, and detector position r are specified. In what fol- 
lows, we restrict our analysis to a semi-infinite medium 
occupying the half-space z > 0, and illuminated by a 
plane monochroma | ti ( j y ave incident at z — 0. For s>l, 
P(r, s) becomes 



P(r,s) 



/ iz 1 



1/2 



exp 



(0) 



Once the field correlation C^(r,r) is known, the 
autocorrelation of intensity fluctuations Csi(r, r) = 
(5I(r,t)SI(r,t + t)) can be found applying the factor- 
ization approximation: Csi(r,r) = \C^(r, t)| 2 £3 

Combining together Eqs. one obtains the nor- 

malized autocorrelation function of multiple-scattered 
wave in a semi-infinite disordered medium: 



9[ L \r,r) 



CV(r, r) 
CV(r,0) 



exp 



a (t) - — 

Lin 



(7) 



where the superscript (L) denotes the linear case, 
a 2 (r) = 3t/(2t ) + £ 2 /L 2 a7 and L a = {££ a /i) 1/2 » £■ 
For the diffusely reflected wave, we assume z ~ £ and get 



g[ L \e,T)^g^>(T)=exp+-a(T) + 



(8) 



From here on, we will use gi(r, t) to denote the nor- 
malized autocorrelation function at a point r inside the 
medium, while <?i(r) — for the normalized autocorrela- 
tion function of diffusely reflected wave. The superscripts 
(L) and (NL) will be used to distinguish between linear 
and nonlinear cases. 



Now we turn to the nonlinear medium. Strictly speak- 
ing, the method of path integrals cannot be applied for 
the analysis of Eq. (|l|), once £2 ^ 0. The failure of the 
path-integral technique follows from the fact that this 
approach relies on the superposition principle which is 
not valid for waves in nonlinear media. However, if the 
nonlinearity is weak [which is ensured by the first two in- 
equalities of Eq. (g)], Eq. (Q) is still approximately valid 
provided that its main ingredients P(r, s) and (A(^ 2 (r)) s 
are computed with account for nonlinear effects. To sim- 
plify such a calculation, we assume that the nonlinearity 
is sufficiently weak to validate the following two assump- 
tions: 

(i) Propagation of waves in a weakly nonlinear disor- 
dered medium is diffusive with a mean free path £ 
unaffected by the nonlinearity. This implies that 
nonlinear refraction is negligible at distances of or- 
der £, and consequently, that An 2 k£ <C 1, where 
An = 712/0 j and Iq is the average intensity in the 
absence of absorption. This assumption is an al- 
ternative formulation of the fact that the role of 
nonlinearity is much less significant than that of 
disorder [see also the first two inequalities of Eq. 

(§]■ 

(m) Intensity of the third harmonic remains always 
much smaller than the intensity of the fundamen- 
tal wave. This implies either that ^(r, t) is con- 
sidered as a complex quantity [in this case, Eq. 
(0) is a nonlinear Schrodinger equation, |^(r, t)| 2 
is time-independent, and the third harmonic is not 
generated at all], or that the medium has a suf- 
ficient degree-pi-idispcrsion for the phase match- 
ing conditiorOEa to be violated: |Afc| £ 1 with 
Ak = &3 — 3k and k% being the wave number at 
frequency Soj. 

Assumption (i) allows us to consider the path distribu- 
tion P(r, s) being unaffected by nonlinearity. The only 
object to be recalculated with account for nonlinear ef- 
fects is then (A</? 2 (t)) s . Before going into an explicit 
calculation of (A</j 2 (t)) s , we devote the next section to 
a brief derivation of some important results for linear 
medium. 



III. CORRELATIONS IN A LINEAR MEDIUM 

As indicated above, we consider a monochromatic 
plane wave incident at the surface z = of a semi- 
infinite medium occupying the half-space z > 0. The 
average intensity at £ can be then found in the diffu- 
sion approximation:!'^ (I(r,t)) = I exp(— z/L a ). The 
spatio-temporal correlation function of the field is given 
by a solution of the Bethe-Salpeter equation (see Ap- 
pendix [A| for details of the calculation) : 
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CV(r,Ar,r) = (^(r - Ar/2, i)V*(r + Ar/2, t + r)) 
sin(fcAr) 



fcAr 



■ exp 



Ar z 



(9) 



Let us now consider the correlation functions of in- 
tensity fluctuations. In addition to SI(r,t), which is 
the deviation of intensity from its average value, it is 
convenient to define A/(r,i, r) = I(r,t + r) — I(r,t), 
which is the change of the local intensity during the time 
interval r. While (SI{r,t)) = (AJ(r,t,r)) = 0, the 
correlation functions Csj(r, Ar, r) and C&j(r, Ar, t) = 
(A/(r - Ar/2, t, r) A/(r + Ar/2, t, r)) for Ar < £ can be 
found in the factorization approximation: 

C 5/ (r,Ar,T) = |CV(r,Ar,r)| 2 , (10) 
C A /(r,Ar,-r) =2[C 5J (r,Ar,0)-C 5J (r,Ar,T)]. (11) 



Both correlation functions (|lfj| ) and (11) decrease ex- 
ponentially with Ar/£, and thus become negligible for 
Ar > £. Intensity correlation persists, however, even for 
two points separated by a distance Ar > £. This corre- 
lation is due to the diffusive nature of wave transport in 
the medium-and can be found either using the Lantmviii 
approachoEil or applying diagrammatic methodsEdO. 
We give the details of calculations in Appendix [B|, the 
final results are 



C 5/ (r,Ar,r) 



(k£) 



-r 

2 2 



dKK 



K, ^K 2 +£ 2 /L 2 ,-,—,a(r) 



xQ 

x J (KAR/e) , 
C A /(r,Ar,r) = 



(12) 



(k£) 



2 2 



x AQ 



K, ^K 2 +£ 2 /L 2 al 



dKK 

o 

z Az 



V i 



a(r) ,a(0) 



x J (KAR/e) . 



(13) 



Here we use the cylindrical coordinates: r = {R, z}, Jq 
is the Bessel function of zeroth order, the function Q is 
defined in Appendix |b| and 



AQ[---,a(r),a(0),} = Q[---,a(0)} 
-Q[---,a(T)} 



(14) 



Due to a rather complicated structure of the function Q 
(see Appendix [b]), further calculations can be done only 
approximately. In the case of z ± Az/2, Ar -C £/a(r), 
we get 



C«(r,Ar, r) 



7-2 
i 



4 (k£) 2 Zg 



1 



x exp 



-2a (r) 



(15) 



where z g is the geometrical average of z-coordinates of 
the two points for which the correlation is computed: 



z g = Jz 2 - Az 2 /4, and x = Ar /( 2 %)- For a(r) = 0, 
Eq. (|lj) is exact. If \ < 1: the correlation behaves essen- 
tially as 1/Ar, while for x > 1 it becomes proportional 
to z 2 /Ar 3 . For the correlation function of Eq. ( |l3|) we 
find: 



C* A /(r,Ar,r) 

(A [a (r)-o(0)] (z/£)G 6I (t, Ar,0), 
~ <^ 2a{T)(z/£) < 1 

l2C«(r,Ar,0), 2a(r)(z/^) > 1 



(16) 



IV. DEPHASING OF WAVES IN A NONLINEAR 
MEDIUM 

Consider a single wave path of length s going from the 
source of waves to some point r. The phase acquired by 
a wave traveling along such a path can be written as 



(f(t) = / fc n[r(si),t] dsi, 



(17) 



where the integration is along the path, and n(r, t) = 
no + ri2l(r,t). The squared difference A<p(t, r) = ip(t + 
t) — (p(t), averaged over various realizations of disorder, 
and over all possible paths of length s, is found directly 
from Eq. (0): 



(18) 



<^ 2 M) S = £< A ^))1 J) > 

3=0 



where the four contributions corresponding to j = 
0, . . . , 3 originate from different physical processes. Be- 
low, we give explicit expressions of these terms and dis- 
cuss their origin. 

,(o) . 



The first term in Eq. (18), (Aip 2 (t)) s is the linear 
term given by Eq. (|^). The next three terms, namely 
the terms corresponding to j — 1, 2, 3, are absent in the 
linear case and only appear because of nonlinear nature of 
wave interaction with the medium. Explicit expressions 
for these terms are 



/2n 2 t 
\ no£ r j 
2n 2 t 
n ro 



(I(r)) dsi 
•s 



((/(r))) s/ - 



(A V 2 (r))f =(— konl f C A i(r,0,T)d Sl 

\ n Jo 



— ko£n 2 2 (C A i(r,0,r)) s - 



(A^(t)) 



(3) 



"o 

k 



o n i [ [ C A /(r, Ar, r)ds 1 ds 2 
Jo Jo 

(M) 2 ^(C A /(r,Ar,r)) s (^) 2 . 



(19) 



(20) 



(21) 



Here the integrations are assumed along wave paths of 
length s ^ £ [in Eq. j2l|), both integrals are along the 
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same path] . Eq. (|2C|) originates from the short-range cor- 
relation of intensity fluctuations [see Eqs. (^)-(ll)]: 



(A^(T)>f = k 2 n 2 2 (C AI (r,0,r)) 

, n 2 



ds / d(As) 
o J-£ 



sin(fcAs) 



kAi 



exp 



A, 



(22) 



where the wave path is assumed to be ballistic at dis- 
tances shorter than I. Eq. (|22l) r educes to (|20f ) for k£ 3> 1. 
Next, the term given by EqT(|2l|) is due to the long-range 
correlations of intensity fluctuations [see Eq. (|l3|)]. The 
averages entering into the right-hand sides of Eqs. (|l9|)- 
(E2) are 



«/(r)» { 



d 3 rp s (r) (7(r)), 



(23) 



(C A /(r,0,r)) s = f d 3 rp a (r)C AI (r,0,T), (24) 

(C AI (r,Ar,r)) s = J d 3 n J d 3 v 2 Ps (r u r 2 ) 

C ai (t,Av,t), (25) 

where the integrations are over the volume of the disor- 
dered medium, and ri, 2 = r ± Ar/2. In Eqs. (p3|)-(p5|), 
p s (r) is the probability density for a path of length s, 
going from the source of waves to the point r, to pass 
through the vicinity of r, and p s (ri,r2) is the probabil- 
ity density for the path to pass consequently through the 
vicinities of ri and r 2 . These two "path distributions" 
should be calculated with account for a particular geom- 
etry of disordered sample and source of waves. Once the 
geometry is fixed, the calculation is straightforward. For 
the case of a plane wave incident upon a semi-infinite dis- 
ordered medium, the calculations of p s (r) and p s (r l7 r 2 ) 
are presented in Appendix ^|. 

Let us discuss briefly the physical origins of nonlinear 
contributions to the dephasing given by Eqs. (|l~9|)-(f2l"|). 

(A<p 2 (r))[, :L ' ) describes the change of the effective wave 
number in the nonlinear medium: k — kgUQ — > £;(r) = 
ko[no + ri2 (7(r, t))]. This contribution can be either pos- 
itive or negative, depending on the sign of n 2 , but its ab- 
solute value is always much smaller than (Aip 2 (t))^ 



as 



long as An <C no. (Af/? 2 ^))^ 1 ' can therefore only cause 
a small correction to the linear correlation function. The 
next contribution, (A</? 2 (t)) s , originates from fluctua- 



in their turn, caused by the scatterer motion. The im- 
portant point is that the scatterer motion is not the only 
possible reason for the fluctuations of intensity with time. 
Weak, spontaneous fluctuations of I(r, t) (due to thermal 
fluctuations of various parameters, vibrations in the ex- 
perimental setup, fluctuations of the incident wave, etc.) 
are inevitable in real physical systems. Eqs. (^0|) and ( pl| ) 
provide a mechanism for this spontaneous and generally 
weak fluctuations to affect the dephasing (A<^ 2 (r)) s and, 
consequently, the temporal correlation function of scat- 
tered wave. 



V. PERTURBATION THEORY 

As stated in the title, the present paper is devoted 
to weakly nonlinear disordered media. We limit our- 
selves to a weak nonlinearity, as otherwise the problem 
becomes too involved. Above, we have already men- 
tioned that we assume An 2 k£ <C 1, and that this con- 
dition allows us to consider the transport of average in- 
tensity to remain unaffected by the nonlinearity [assump- 
tion (i) of Sec. || . This allows us to use "linear" results, 
(I(r)) = I Q exp(-z/L a ) and Eqs. (Jc|), Q of Appendix 
g, for (J(r)), Ps (r), and p s (n,r 2 ) in Eqs. 
It seems then natural to assume that C A i(r, 0,r) and 
C A i(r, Ar, r) are also close to their linear values. We 
can therefore replace these correlation functions in Eqs. 
( p4| ) and ( p5| ) by the exp ressions found in Sec. III. Then, 
making use of Eqs. (|op-(|l|), and performing necessary 
integrations, we obtain from Eqs. (|l9|)-(21 ): 



<A^ 2 (r)>; 13 = 2An-\l-H 



TO 



(A( /3 2 (r))f ) = 2Trk eAn 2 jff 



- H 




(Atp 2 (T))f ] = 6An 2 S a(r)V^,a(0)V^ 

3/2 




(26) 



(27) 



(28) 



Here H(x) — ^/ttx exp(x 2 )[l — Erf(x)] and no = 1 is 
assumed for simplicity. The function S(u, v) in Eq. d2§| ) 
is 



tions of the local intensity, while (A<^ 2 (r))^ is due to S{u,v) = 9 / dRR / dK K / dz d(Az) 

the lonc-ra.up'e correlation of these fluctuations. An im- Jo JO JO JO 



the long-range correlation of these fluctuations. An im 
portant difference between the linear term (^|), the first 
nonlinear term (JTsj) , and the two last nonlinear terms 
(f^of), (|2l]), is that the latter terms do not depend ex- 
plicitly on To. The terms given by Eqs. (^0|) and ( |2~l| ) 
are determined by the intensity fluctuations, and not by 
the scatterer displacements. This might seem to be a 
meaningless statement, as the intensity fluctuations are, 



f(z,Az,R) 



x f(z, Az, R)AQ(K, VK 2 +v 2 , z, Az, u, v) 
x Jo(KR), (29) 
2z + VR 2 + Az 2 



VR 2 + Az 2 



x exp 



-- (2z+y/B? + Az' 



G 



2z + y/R 2 + 4z 2 



x exp 



Vi? 2 + 4z 2 
3 



2z + \]R 2 + Az 2 



(30) 



Unfortunately, integrations in Eq. (|29|) cannot be per- 
formed in the general case. We find, however, the follow- 
ing approximate results: 



S(u, v) 




v )v 



-3 



u — V < 1, V < 1, 

it — u > 1, v < 1, 
W — v < 1,« > 1, 

u — V > 1, V > 1. 



(31) 



Here numerical factors of order unity are omitted before 
each of the four asymptotic expressions. 

While approximate, the above results enable one to 
compute the temporal autocorrelation function g/^ W) 
of diffusely reflected wave numerically using Eqs. (|4|), (g), 
and (|f): 



gf L \r) = W[{A l p 2 (r)) s } /W[0], 
W[(A^(r)) s ] 



P{£, s) exp 



-§<V(t)>, 



ds. 



(32) 
(33) 



We present the results of the calculation in Fig. [l] for 
fixed An, kg£, and the values of the inverse absorption 
length £/L a indicated near each curve. The "linear" 
correlation functions, corresponding to the same three 
values of the absorption length, and to An = 0, are 
shown by dashed lines. For weak absorption (£/L a = 0, 
10~ 3 ), the initial (i.e. at short correlation delay times 
r) decrease of the "nonlinear" autocorrelation function 
is much faster than that of the linear one. Hence, our 
perturbation approach fails at short times. Indeed, the 
results of this section are based on the assumption that 
Ca/(i",0, t) and Ca/(i", Ar, r) in the nonlinear medium 
are close to their values in the linear one. Applying 
the factorization approximation, we find that this im- 
plies that 1 - \g[ NL) {T)\ 2 « 1 - |5i L) (r)| 2 - Although 
this condition seems to hold well for sufficiently large 
r, it can be violated at small r, where, as follows from 
Fig. [D, 1 |Si ( r )| 2 can become much larger than 
1 — Iffi (t)| 2 , if absorption is sufficiently weak. 

To estimate the region of validity of our perturbation 
approach, we require that the linear contribution to the 
dephasing (A<p 2 (r)) s given by Eq. (||), should be consid- 
erably greater than the sum of nonlinear contributions 
[Eqs. (Pq)-(^)]. As the longest path length contributing 
to the integral of Eqs. (Q) and ([53]) is s ~ £/a 2 , we ob- 
tain the conditions of validity of the perturbation theory 
in the form: 



An 2 a(rT 2 [k £ + a(r)" 1 ] < 1. 



(34) 
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FIG. 1. Normalized temporal autocorrelation function of 
a wave diffusely reflected from a semi-infinite nonlinear dis- 
ordered medium, calculated using perturbation theory for 
An = 10~ 4 , ko£ — 100, and the values of £/L a indicated 
near each curve (solid lines). Dashed lines show correspond- 
ing results for a linear medium (An = 0). 

Since ol{t) is an increasing function of its argument, 
and a(0) = i/L a , condition An 2 (L a /£) 2 [k £ + Lji] < 1 
ensures Eq. ( |34| ) at any r. It is the case for the up- 
per curve of Fig. |l|, corresponding to £/L a — 5 x 10~ 3 . 
For such a strong absorption, the perturbation theory 
is valid at any r, and the nonlinear autocorrelation 
function is close to the linear result. By contrast, if 
An 2 (L a /£) 2 [k £ + L a /£] > 1, the perturbation approach 
can be applied only for sufficiently long correlation times 
t > t c , where the critical time r c is determined by Eq. 
(|34|). For the lower curve of Fig. [j], corresponding to 
£/L a = 0, we find (t c /t ) 1/2 « 2 x 10~ 3 . It is worth- 
while to note that in the absence of absorption, even an 
infinitely small nonlincarity (An — > 0) suffices for the 
perturbation theory to fail at t/tq — > 0. Our condition 
of validity of the perturbation approach (H) is consis- 
tent with the result of Spivak and ZyuzinJHI who have 
shown that the perturbation analysis of the sensitivity 
of speckle pattern in a nonlinear disordered medium to 
changes of scattering potential fails for An 2 (L/£) 3 > 1, 
where L 3> £ is the typical size of the medium. 



VI. SELF-CONSISTENT ANALYSIS 

As demonstrated in the previous section, the pertur- 
bation theory fails to describe the temporal autocorrela- 
tion function of wave diffusely reflected from a nonlinear 
medium for r < r c , where t c is defined by Eq. (|34|). To 

calculate g[ NL ^ (r) at such short times, one has to use a 
non-perturbative, self-consistent analysis. However, the 
possibility of performing such an analysis is considerably 
limited by the mathematical complexity of the considered 
problem. To make the self-consistent analysis possible, 
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we adopt the following two additional assumptions: 

(Hi) The statistics of a wave field scattered in a 
weakly nonlinear disordered medium, is close to 
Gaussian. Consequently, the factorization ap- 
proximation holds in a weakly nonlinear medium: 
C si (t,t)~\C^(v,t)\ 2 . 

(iv) The functional form of the r-dependence of 
g[ NL \r, t) is the same as that of g[ L \r,r): 

g[ NL \r,r) = exp[— (3(t) z/£], where (3(t) is some 
unknown function, which can depend not only on 
r but also on other parameters of the problem 
(namely, on £/ L a , k§£, An). 

Strictly speaking, the above assumptions define a sort 
of perturbation theory, but now we do not limit the values 
of deviations of intensity and field correlation functions 
from their values in the linear case. Instead, we assume 
that the nonlinearity does not cause significant modifi- 
cations of the statistics of scattered waves [assumption 
(Hi)] and of the functional form of the field correlation 
function [assumption {iv)]. Note that now the autocorre- 
lation function g[ NL ^ (r, r) can deviate significantly from 

ffi (r, t). Condition (iv) fixes the functional form of 
this deviation, but implies no constraints on its absolute 
value. 

Obviously, both the assumptions (Hi) and (iv) require 
the nonlinearity to be weak. Assumption (Hi) is justi- 
fied under the same conditions as (i) (see Sec. [0]), since 
the Gaussian statistics of the total scattered wave field 
?/>(r, t) is a consequence of the complete randomization 
of phases of partial waves arriving at r. The reason for 
the randomization is that the typical distance £ between 
individual scattering events in a multiple-scattering s 
quence is much larger than the wavelength (k£ 3> 1) 
Obviously, such a mechanism of phase randomization is 
equally effective in both linear and weakly nonlinear me- 
dia, as long as (i) holds. 

To justify the ansatz of assumption (iv), we apply Eq. 



and write jj ' as 



(NL), s 

Si (r,T) 



1/2 



exp 



„3/2 



exp 



I ds, 

Ms ' 



(35) 



where we assumed that s/£ a + (1/2) (Aip 2 (r)) s increases 
monotonically with s and becomes of order unity at 
s = Si(t) 3> so = z 2 /(2£). After performing in- 
tegration, Eq. ( p5| ) can be approximately rewritten as 
exp(z/L a ){l - z~$/(ir£s l )] 1 i 2 } ~ exp[-j8(r) z/£] with 
f3z/£ ~ z(£s\)" x l 2 — z/L a <C 1. In the opposite limit of 
/3z/£ 3> 1, <7i (r,r) vanishes and the functional form 
of its z-dependence is of no importance. Anyway, it will 
be seen from the following that the exact functional form 



of the r-dependence of <7i(r, r) is not of crucial impor- 
tance, since gi is integrated over the whole medium dur- 
ing the calculation of the correlation function of diffusely 
reflected wave. 

Making use of Eqs. (|l9|)-(|2"l|) and relying on the as- 
sumption (Hi), (iv), we recalculate the nonlinear contri- 
butions to the dephasing (A^ 2 (r)) s . Due to the assump- 
tion (iv), the result can be obtained from Eqs. j26|)-(|2^) 
by a simple substitution: a(r) — > j3(r) +a(Q). Recalling 
that a(Q) = £/L a , we obtain: 



(A^(r)). W =2A 



(Aip 2 (T))f } = 2nk £An 



- H 




(A^ 2 (r)) 



(3) 



X < 



£ 

6An 2 

r & x ( s /£) 2 , 

0y/tj£<l,(l/L a ) 
(s/£f/ 2 , 

Pyftji>\,{l/L a ) 

(3 x (L a /£) 3 x yfi/t, 

Vyfifl < 1, (£/ L a ) 
(La/£) 3 , 

/3y/I/l>l,(£/L a ) 



(36) 



(37) 




(38) 



As follows from Eqs. (p6|)-(38), the nonlinear dephas- 
ing depends on the unknown function [3. Since the tem- 
poral autocorrelation function of diffusely reflected wave 

ffi ( T ) = ex P[ — P{ T )] is determined by (3 as well, Eqs. 
( p2| ) and (|33|) allow us to formulate a self-consistent equa- 
tion for (3: 



cx P [-/?(t)]=F[/3(t)], 
where F(J3) = W [(A<^ 2 (t))J /W[Q], W 



(39) 
is defined 

by Eq. ([33), and (A(/3 2 (t)\ is a sum of terms given by 
Eqs. (§), (3^-@. Eq. @ is the main result of our 
self-consistent analysis. Although the functional form 
of F(j3) is rather complicated, and Eq. (|39| ) cannot be 
solved analytically, the numerical solution is straightfor- 
ward and can be carried out for any values of An, £/L a , 
kg£ and r. Eq. (|3^) can be considered as a self-consistent 
equation for the autocorrelation function of diffusely re- 
flected wave, since (3(t) and g[ NL \r) are directly related. 
It is worthwhile to note that in the absence of nonlinear- 
ity (An = 0), Eq. @ yields (3(t) = a(r) - £/L a , and 
hence Eq. (|J) is recovered for g[ L ^ (t) . 



VII. RESULTS AND DISCUSSION 



We start the analysis of Eq. ( |39[ ) from the case of im- 
mobile scatterers, taking a limit of t/tq — ► 0. We de- 
note the autocorrelation functions corresponding to this 
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limit by g[ L ^^ NL \o + ) in order to distinguish them from 
g[ L ^'^ NL \o), which correspond to t = 0@ Obviously, 

g\ L \o + ) = fli (0) = 1- In a nonlinear medium, Eqs. 
(B71) and (Bq) can still contribute to the dephasing even 
for t/tq — > 0. These contributions are insensitive to the 
sign of An. A corresponding value of /3(0 + ) and, conse- 
quently, of the field autocorrelation function g[ NL \o~D 
can be found by solving Eq. ( ^9| ) numerically. In Fig ,J2, 
we plot the left-hand and the right-hand sides of Eq. (B9) 
for fixed |An| = 1CP 4 , k^l = 100, and for several values 
of £/L a . If absorption is strong (l/L a > 3 x 10~ 3 for 
considered |An| and ko£), Eq. ( p9| ) has a unique solution 

/3(0 + ) = which corresponds to g[ NL \o + ) = 1- How- 
ever, a second solution /3(0 + ) > appears for sufficiently 
weak absorption (£/L a < 3 x 10~ 3 ). The point of ap- 
pearance of the second solution is a bifurcation point of 
Eq. (|39|). To choose the solution realizable in a physical 
system, we note that the first solution [/?(0 + ) = 0] exists 
only for t/tq = 0, and disappears at finite t/tq, since 
F(0) < 1 for t/t > 0. This solution is therefore in- 
accessible by continuity, and "unstable" with respect to 
small scatterer displacements. The physically realizable 
solution should represent the limit of /3(r) for t/tq — > 0, 
which is given by the second solution of Eq. (p^). It is 
therefore this solution which one expects to be realized 
in a real physical system. 

The fact that Eq. ([39]) can have a positive solution 
/3(0 + ) > is a very important issue, since /3(0 + ) > 

leads to g[ NL) {0+) = exp[-/3(0+)] < 1. A value of 
the temporal autocorrelation function which is less than 
unity is commonly associated with temporal fluctuations 
of scattered waves. In the considered case, however, the 
reason for these fluctuations is not the motion of scat- 
terers, as the limit of t/t — > corresponds to immo- 
bile scatterers. The fluctuations are spontaneous and 
represent a clear signaturq-af instability of the multiple- 
scattering speckle pattern.E9 

Despite a rather complicated structure of the function 
F((3) in Eq. (|39|), a relation between the parameters of 
the problem corresponding to the onset of the speckle 
pattern instability can be found analytically. As follows 
from Fig. |, the initial [at /3(0+) = 0] decay of F[/?(0+)] 
should be faster than the decay of exp[— /3(0 + )], for the 
second solution of Eq. fl39| ) to appear. A surface in a 
three-dimensional space of the problem parameters An, 
i/L a , ko£, separating the stable [/?(0 + ) = 0] and unstable 
[/3(0 + ) > 0] regions, is therefore given by the equation 
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FIG. 2. Graphical solution of Eq. (B9|) at r/r -> 0. Solid 
lines show F\J3(0 + )] for \An\ = 10,k £ = 100, and the 
values of £/L a indicated near each curve. Dashed line is 
cxp[— /3(0 + )]. If absorption is weak {i/L a = 0, 10" 3 , and 
2xlO _3 ),Eq. @ has two solutions /3(0+) =0and/3(0 + ) > 0, 
while for strong absorption (£/L a = 3x 10 -3 and 5 x 10 -3 ), 
the second solution disappears. 

Recalling that L a /£ ^> 1, k £ ^> 1 is assumed, we 
obtain from Eq. (40): 



p = An ( -J 



k £+^ 



1. 



(41) 



where we introduce a control parameter p, and a numeri- 
cal factor of order unity is omitted. If p < 1, the multiple- 
scattering speckle pattern is stable [g[ NL \o + ) = 1], while 
for p > 1, an instability shows up leading to g[ NL \o + ) < 
1. A striking feature of Eq. ( pl| ) is that p can become 
larger than unity even for very small |An|, provided that 
the extensive parameter L a /t is large enough. Our con- 
dition of the speckle pattern stability p < 1 agrees with 
the condition of validity of the perturbation theory de- 
veloped in Sec. |y| [Eq. (gj)], evaluated at t/tq = 0. This 
readily explains the failure of the perturbation theory for 
short correlation times and weak absorption: perturba- 
tion approach is not suitable for description of unstable 
regimes. Moreover, the fact that the condition of valid- 
ity of the perturbation theory and the speckle pattern 
stability condition agrees indicates that the additional 
assumptions (m) and (iv) of Sec. VI are not essential for 
obtaining unstable regimes. 
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FIG. 3. "Bifurcation diagram" for a wave scattered in a 
semi-infinite nonlinear disordered medium. Solid lines show 
the temporal autocorrelation function of diffusely reflected 
wave at t/to — > (immobile scatterers) for ko£ = 100 and 
the values of i/L a indicated near each curve (dashed line is 
for £/L a = 0). The threshold values of |An| following from 
Eq. ( |4ll ) are shown by arrows. g[ NL \o + ) < 1 corresponds to 
spontaneous fluctuations of scattered wave, which is a mani- 
festation of the speckle pattern instability. Dotted line is the 
linear result g[ N (0 + ) = 1. 

To illustrate our self-consistent theoretical framework, 
we solve Eq. (39) numerically for t/to — * 0, and plot the 
resulting temporal autocorrelation function of diffusely 
reflected wave, g[ NL \o + ), in Fig. [| As discussed above, 

<?i (0 + ) < 1 corresponds to spontaneous fluctuations 
of scattered waves which is a manifestation of the speckle 
pattern instability. It follows from Fig. 0, that in accor- 
dance with Eq. (41), an infinitely small |An| is sufficient 
to make the speckle pattern unstable in the absence of 
absorption, while a certain threshold degree of nonlin- 
earity is required to destabilize the speckle pattern in a 
dissipative medium. In the absence of absorption, Eq. 
( |39l ) always has two solutions: /3(0 + ) = and /3(0 + ) > 0, 
corresponding to g[ NL \o + ) = 1 and g[ NL \o + ) < 1, re- 
spectively. As discussed above, it is the second solution, 
shown by a dashed line in Fig. ||, which is the physically 
realizable one. In contrast, if i/L a ^ 0, |An| should be 
greater than some threshold value for the second solution 
/3(0 + ) > 1 to appear. Threshold values of |An| following 
from Eq. (41) are shown in Fig. |3|by arrows. 

For a nonabsorbing, clastically scattering medium 
(l/L a = 0), the value of g[ NL} (0+) can be estimated ana- 
lytically. Indeed, at t/tq — > the principal contribution 

by <A^ 2 (r))f for s/i < (k i) 2 , 
for s/i > {k £) 2 . This allows us 



to (A(^ 2 (r)) is given 
and by (A^ 2 (r))[ 3) 



to put (A(p 2 (t))^ 



(A^ 2 (t)) 

,( 3 ) 



(2) 



if s/i < (k Q i) 2 , and 



{A<p 2 (t)) s a (Av? 2 (t)); ; if s/i > {k i) 2 . Integration in 



Eq. ( |33D can be then carried out, and Eq. ( p9| ) is easily 
solved, yielding 



1 



2\An 



2/3 



3|An| Al, 



|An| < (k Q i)~ 3 / 2 , 
\An\ > {k l)- 3 / 2 . 



(42) 



This result agrees well with the numerical solution of Eq. 
([39|) presented in Fig ||. 

The physical origin of the instability of speckle pat- 
tern for p > 1 can be revealed by realizing that the 
system "coherent wave + nonlinear disordered medium" 
has a positive three-dimensional feedback. In a nonlin- 
ear medium, an infinitely small perturbation of the wave 
intensity /(r, t) produces a change of the local refractive 
index, which alters the phases of waves propagating in 
the medium and, consequently, affects their mutual in- 
terference. Since it is this interference which determines 
I(r,t), the loop of the feedback is closed. For p 1, 
the feedback is sufficiently strong to compensate for the 
(diffusive on average) spreading of the initial intensity 
perturbation, and the speckle pattern /(r, t) is unstable. 
It is worthwhile to note that unstable regimes are not 
exceptional in nonlinear wave systems and, in particular, 
in optical systems (see, e.g., Refs. f72|-f76|). 
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FIG. 4. The same as Fig. ^, but using a self-consistent ap- 
proach instead of the perturbation theory. For the two lower 



curves, g[ NL \0 + ) < 1 and the speckle pattern is unstable. 



Our Eq. (|39|) is in no way limited to the case of 
t/tq — > 0, and can be used to compute the tempo- 
ral autocorrelation function of diffusely reflected wave 
at any t/tq > 0. In the latter case, one should take 
into account all the four dephasing terms given by Eqs. 

(|36|) ( p8|) . The results of the numerical solution 
of Eq. (]39|) are shown in Fig. || For weak absorption 
(i/L a = 0, 10~ 3 ), p > 1 and the speckle pattern is un- 
stable [g[ NL) {0 + ) < 1]. As r/r increases, the difference 
between the "nonlinear" and "linear" curves becomes less 
pronounced. For strong absorption (i/L a = 5 x 10 -3 ), p 
becomes less than unity and stability of the speckle pat- 
tern is recovered [g[ NL \o + ) — !]■ We remind that the 
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temporal autocorrelation function g[ (t), correspond- 
ing to a linear medium, always equals 1 for t/t = 0, as 
shown by dashed lines in Fig. [|. 
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FIG. 5. Comparison of the normalized temporal autocorre- 
lation functions calculated using a self-consistent (solid lines) 
and perturbation (dashed lines) theories. Dotted vertical line 
indicates the limit r c of validity of the perturbation theory for 
l/L a = 0. The perturbation theory is valid only for r > r c , 
while for r < r c the perturbation result deviates significantly 
from the self-consistent solution. For sufficiently strong ab- 
sorption (l/L a = 5 x 10~ 3 ), the perturbation theory holds 
at any r: the perturbation and the self-consistent results are 
indistinguishable . 



(b). The onset of the speckle pattern instability occurs 
when the control parameter 



p = An 2 



2 r 



hot- 



La 

i 



(43) 



becomes of order or larger than unity. The speckle 
pattern is stable for p < 1. 

The instability of the multiple-scattering speckle pattern 
manifests itself in spontaneous fluctuations of the scat- 
tered wave field and intensity. The following features 
are characteristic for the development of the instability. 
First, the development of the instability is independent of 
the sign of nonlinearity. This is not common for nonlin- 
ear waves since the instability is often due to_sehjfocusing 
phenomena, which only occur for n-i > 0.L51Z3 The in- 
stability of waves in a disordered nonlinear medium has 
nothing to do with the self-focusing, and occurs at rel- 
atively weak nonlinearities, when the self-focusing can 
be neglected. Second, in the absence of absorption, the 
speckle pattern is unstable for any (even infinitely small) 
value of the nonlinearity strength |An|, while in a dis- 
sipative medium |An| should exceed a certain threshold 
value for the instability to show up. Finally, the insta- 
bility results in a value of the autocorrelation function 
g[ NL \T) of scattered wave, which is smaller than 1 for 
t/tq — > (i.e., in the absence of scatterer motion). This 
is a clear signature of spontaneous fluctuations of the 
multiple-scattered speckle pattern, and should be observ- 
able in experiments. 



It is instructive to compare the results obtained using 
the perturbation theory of Sec. ^ and the self-consistent 
approach of Sec. VI. Such a comparison is shown in 
Fig. ||. The two upper curves, corresponding to rela- 
tively strong absorption (£/L a — 5 x 10~ 3 , p < 1), are 
almost indistinguishable, which means that for p < 1, 
the perturbation theory works very well. In contrast, 
for p > 1 (see the two lower curves of Fig. |5| correspond- 
ing to £/L a = 0), the perturbation and the self-consistent 
curves are close only at the right from the dotted vertical 
line, showing the minimum time r c at which the pertur- 
bation theory is valid [see Eq. (Q)]. For r < r c , the per- 
turbation and the self-consistent results disagree signif- 
icantly, which confirms our conclusion about the failure 
of the perturbation approach at short correlation delay 
times. 



VIII. CONCLUSION 

We are now in a position to answer the two central 
questions formulated in the introductory section: 

(a). The phenomenon of multiple scattering is capable 
of providing a positive feedback for a coherent wave 
propagating in a nonlinear disordered medium. 
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APPENDIX A: SPATIO-TEMPORAL FILED 
CORRELATION IN A LINEAR MEDIUM 

In this appendix, we provide a derivation of the 
spatio-temporal correlation function, C^,(r, Ar, t) = 
(ip(r- Ar/2,£)V*(r + Ar/2, i + r)), of a random field 
ijj{v,t) in the bulk of disordered medium. Starting from 
the linear wave equation [Eq. ([!]) with e2r-^r£lj|r.we obtain 
the Bethe-Salpeter equation in the formOSB 

CV(r,Ar,r) = (V>(r - Ar/2, t)) (V*(r + Ar/2, t)) 

dr a / dr b / dr c / dr d G{v - Ar/2, r„) 



xG (r + Ar/2,r b )[/( 

*"a: Tfr, r c , r^) 
x C^[(r c + r d )/2,r d - r c ,r], 



(Al) 



where the integrations are over the volume of disordered 
medium, G is the average Green function of the linear 
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wave equation, and U is the irreducible four-point vertex. 
Far enough from the medium boundaries, we can replace 
G by its value in the infinite medium: 



G(ri,r 2 ) 



1 



Air ri - r 2 



■ exp 



ik 

2£ 



l r i - r 2 | 



(A2) 



Now we assume that point scatterers in the medium 
undergo Brownian motion, and that the correlation func- 
tion of the dielectric function fluctuations is given by Eq. 
(||). For t < To, we can neglect the r-dependence of U 
in Eq. (Al), which in the limit of k£ 3> 1 becomes 



U(r a ,r b ,r c ,r d ) 

= -f S i r a ~ r c)S{r b - r d )S(r a - r b ). 



(A3) 



Eq. (Al) then reduces to 
CV(r,Ar, r) 



4tt 

T 



dv a G{v- Ar/2,r ) 



xG (r + Ar/2,r a )CV(r a) 0,T), 



(A4) 



where we assumed the coherent field (ip(r,t)) to be neg- 
ligible. 

Now we remind that C^,(r a , 0, r) = C^(r , r) is a slow- 
varying function of r a . We can therefore pull it out of 
the integral in Eq. (A4), taking its value at r. Performing 
the remaining integral, we find 



CV(r,Ar,r) 



sin(fcAr) 
kAr 



■ exp 



This equation relates the spatio-temporal correlation of 
the field with its purely temporal correlation. Eq. (A5) 
holds for any sample geometry and source distribution, 
far enough from boundaries and sources, and for r <C tq. 
In the case of a plane wave incident at the surface z = of 
a semi-infinite disordered medium, occupying the z > 
half-space, 



C/,(r, r) = I exp [-ol{t){z/1)\ 



(A6) 



with a 2 (r) = 3r/(2r ) + (£/L a ) 2 , leading to Eq. (§) i n 
the main text. In the factorization approximation (|10[), 
the square of Eq. ( A5) gives the short-range correlation 
function of intensity fluctuations. 



APPENDIX B: LONG-RANGE 
SPATIO-TEMPORAL CORRELATION OF 
INTENSITY FLUCTUATIONS IN A LINEAR 
MEDIUM 

To calculate the spatio-temporal correlation function 
of intensity fluctuations, 

C 5/ (r,Ar,r) = {51 (r - Ar/2, t)6I(r + Ar/JLt-J^)) , for 
Ar > I we generalize the Langevin approach,0OE3 which 



has been initially developed for calculation of purely spa- 
tial correlations. Wc assume that the spatio-temporal 
correlation function of Langevin random sources is de- 
termined by the short-range correlation function of in- 
tensity fluctuations. In the factorization approximation, 
the latter is given by t he square of the field-field correla- 
tion function [Eq. ( A5)]. Hence, the generalized Langevin 
equations read 



D p [V 2 - l/L 2 a ] SI(v,t) = div.UM), 
£(r-Ar/2,i)iW(r + Ar/2,i + r) 



-h — 

- 3 d "» k 2 



|C^(r,r)r«(Ar), 



(Bl) 



(B2) 



where D p = c£/3 is the diffusion coefficient, and c is 
the speed of wave in the medium. In the case of a 
plane wave incident upon a boundary z — of a semi- 
infinite disordered medium, it is convenient to make a 



two-dimensional Fourier transform of Eq. (Bl) in the 
{x, y} plane. Then the equations corresponding to 
(5/(Ki, z\, ti) and <5/(K2, z 2 , i 2 ) are multiplied, and the 
product is ensemble averaged. Further transformations, 
which are equivalent to those discussed in Ref. [0], yield 

(6I(K 1 ,z 1 ,t 1 )6I(K 2 ,z 2 ,t 2 )) = I 2 ^(Ki - K 2 ) 



dz' 



(K 1 -K 2 )G(p,z 1 ,z')G(p, z 2 ,z') 



d d 
+ -^G(p,zi,z)—G(p,z 2 ,z) 
az oz 

where p 2 = K 2 + 1/L 2 , and 



exp 



-2«(r)j 



(B3) 



G(p, zi,z') 



1 

P 



sinh(p x min{2i,2;'}) 



x exp(— p x maxlz!,^'}) 



(B4) 



is the Green function of Eq. (Bl). Evaluating Eq 
and transforming the result back to the real space, after 
lengthy but straightforward algebra we obtain 



C 5J (r,Ar,r) 



(H) 2 ° Ja 



dKK 



x Q 



xJ [K 



.AR 



(B5) 



where 



Q(K,p, C, AC, a) = exp (2p£) x 



2a 2 + K 2 



4a 3 — Aap 2 



[{a 2 - p 2 )(-K 2 + p 2 ) + a 2 {K 2 + p 2 ) cosh(2pCi) 
ap(K 2 +p 2 )smh(2Ki)] / [4aexp(2aCi)(a 2 - P 2 )p 2 ]] 
K 2 \ sinh(p(i) sinh(pC2) 



1 



P 



2exp[2(a+p)C2] (a + p) 
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exp(-K2)sinh(pCi) 
Aap 2 (a + p) 

x {p(p 2 + 2ap - K 2 ) [cosh(K2) exp[-(2a + p)£ 2 ] 

- cosh(pCi)exp[-(2a + p)Ci]] 

+ (p 3 - 2aK 2 - pK 2 ) [sinh(K 2 ) exp[-(2a + p)C 2 ] 

-smh(pCi)exp[-(2a + p)Ci]]}, (B6) 

where Ci = C ~ AC/2 and (2 = C + AC/2. 

Equations (B5), ([B6| ) cannot be evaluated in a general 
form. In the limits of a(T)£,a(r)AR/£ < 1 we can, 
however, approximately replace Q in Eq. ( |B5| ) by 

Q(K,p, £ AC, a) « ^ [exp(-XAC) - exp(-2XC)] 



x exp[— 2a(r)C] 



(B7) 



which yields Eq. (15) in the main text 



APPENDIX C: DENSITY OF WAVE PATHS IN A 
DISORDERED MEDIUM 

Propagation of waves in a disordered medium can be 
interpreted in terms of partial waves traveling along var- 
ious paths inside the medium. The spatial distribution 
of such paths and their relative weights depend on the 
scattering properties of the medium, and on the geom- 
etry of the sample. In the case of multiple scattering, 
the simplest and, at the same time, sufficiently accurate 
model of wave propagation is the diffusion model. Ac- 
cording to this model, wave paths in the medium coincide 
with trajectories of Brownian particles. The probability 
G(r 4 , r 2 , s) for a path of length s to pass from ri to r 2 is 
then given by a solution of the diffusion equation, which 
in the absence of absorption reads:E£Eil 

^G(n, r 2 , s) - ^V 2 G(n, r 2 , s) = *(r x - r 2 )6(s), (CI) 

where I is the mean fre e path. Commonly used boundary 
conditions for Eq. (CI) consist in putting G = at open 
boundaries and V n G = at reflecting boundaries of the 
sample (where V„ denotes the normal derivative of G). 
G(ri, r 2 , s) is called the Green function, or the propaga- 
tor. For a semi-infinite medium occupying the half-space 
z > 0, one finds 



G(ri,r 2 ,s) 



4n£s 



3/2 



exp 



3 
lYs 



--—(AR 2 + Az 2 ) 



exp 



■jf;(A* + *) 



(02) 



{Rj , Zi\ , 



where cylindrical coordinates are used: 
AR = Ri — R 2 , Az = z\ — z 2 , and Z = z\ + z 2 . 

Following Ref. we introduce p s (ri,r 2 ,r 3 ), the den- 
sity distribution of paths of length s, as a number of visits 
of a given site r 2 inside <i 3 r 2 in the ensemble of paths of 



length s starting at ri and ending at r3, over the total 
length of the ensemble distinct paths: 



/O s (ri,r 2 ,r 3 ) 



1 



/ dpG(ri,r 2 ,p) 
Jo 



sG(ri,r 3 ,s) 
x G(r 2 ,r 3 ,s-p). (C3) 



p s (ri, r 2 , r 3 ) describes the probability density for a path 
of a given length s, starting at rj and ending at r 3 , to 
pass through r 2 . This quantity is normalized: 



J d 3 r 2 p s (ri,r 2 ,r 3 ) = 1, 



(C4) 



where the integration is performed over the volume of 
disordered medium. 

As the Green function G is known [Eq. (|C2[)1, the cal- 
culation of p s (ri,r 2 ,r 3 ) is straightforward. For diffusely 
reflected paths, assuming that the first and the last scat- 
tering events take place at z — £, we obtain 



Ps {£,v,l) = p s {v) 



1 6z 
~ATs 



exp 



" is 



(C5) 



where A — > 00 is the surface of the semi-infinite medium. 



Equation (C5) defines the probability density for a dif- 
fusely reflected path of length s to pass through a vicinity 
of some point r = {x, y, z}. 

Generalizing definition of p s , we define the probability 
density for a path of length s starting at ri and ending 
at Y4, to pass consequently through r 2 and r 3 : 



P s (ri,r 2 ,r 3 ,r 4 ) 



dp 



s 2 G(ri,r 4 ,s) J Q Jo 
G(ri,r 2 ,p)G(r 2 ,r 3 ,g) 
x G(r 3 ,r 4 ,s -p- q). 



dq 



The normalization of Eq. (C6) is 



J d 3 r 2 J d 3 r 3/ 9 s (ri,r 2 ,r 3 ,r 4 ) = 1. 
For a semi-infinite medium, we get 
p s (£,r,r',£)=p s (r,r') ' 



(C6) 



(C7) 



A 2 2ir£ 2 s 2 



Z + VAR 2 + Az 2 



VAi? 2 + Az 2 

'' (z + AR 2 + Az 2 



Ms 



x exp 

z + VaiFTz 2 



x exp 



VAR 2 + Z 2 

■Ws{ Z+ ^ 



z 2 



(C8) 



where r = {R, z}, AR = R-R, Az = z-z', Z = z + z'. 
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Although Eqs. (C5) and ( |Cq ) have been found for a 
nonabsorbing medium, it is easy to show that these re- 
sults hold in the presence of spatially-homogeneous ab- 
sorption as well. This stems from the fact that the atten- 
uation of wave in a homogeneously absorbing medium de- 
pends only on the path length, while being independent of 
the p ath geometry. As a consequence, the Green function 
[Eq. (|C^ )1 should be multiplied by a factor exp(— s/£ a ), 
where £ a is the absorption length. This factor, however, 
disappears after the substitution of the Green function 
( |C2| ) in Eqs. ( |C3| ) and (jC6|). Consequently, p s (r) and 
p s (r, r') are independent of l a and remain unchanged. 
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